Bose-Einstein droplet in free space 
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We show that a droplet of a Bose-Einstein condensate can be dynamically stabilized in free space 
by rapid oscillations of interatomic interactions between attractive and repulsive through, e.g., the 
Feshbach resonance. Energy dissipation, which is present in realistic situations, is found to play a 
crucial role to suppress dynamical instabilities inherent in nonlinear nonequilibrium systems. 
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I. INTRODUCTION 

Matter-wave bright solitons — stable solitary waves 
whose density is greater than the background one — 
have recently been realized in a quasi one-dimensional 
(ID) Bose-Einstein condensate (BEC) In ID, the 

quantum kinetic pressure counterbalances an attractive 
interatomic interaction, allowing a stable bright soliton to 
be formed. However, in 2D or higher dimensions, bright 
solitons are always unstable against collapse or expan- 
sion as long as the attractive interaction is constant in 
time 01 ■ 

It has been demonstrated by AbduUaev et al. ^ and 
by us that by rapid oscillations of interaction between 
attractive and repulsive through, e.g., the Feshbach reso- 
nance |[5j, 1^] , a BEC can be stabilized in a quasi 2D space 
without a radial trapping potential. The mechanism of 
the stabilization is similar to that of an inverted pen- 
dulum, in which a bob is stabilized above the vertically 
oscillating pivot 0, 111- However, the stabilization of a 
BEC in 3D free space has been elusive due to complex 
dynamical instabilities inherent in nonlinear nonequilib- 
rium systems. 

In the present paper, we show that a droplet 9] of a 
BEC can be stabilized in 3D free space by using oscil- 
lating interactions. In order to simulate experimental 
situations at finite temperature, we take into account 
the effect of dissipation caused by the interaction be- 
tween a BEC and a thermal cloud. This can be incorpo- 
rated in mean field theory by a non-unitary term added 
to the Gross-Pitaevskii (CP) equation. This model has 
successfully been used to describe damping of collective 
modes (1^ and vortex nucleation ^| . We find that the 
dissipation plays a crucial role in the stabilization of a 
BEC droplet by suppressing the dynamical instabilities, 
and show that the droplet can indeed be stabilized for a 
realistic amount of dissipation. 

This paper is organized as follows. Section ^ intro- 
duces the GP equation with oscillating interactions, and 
shows that the dependence of the dynamics on the oscilla- 
tion frequency can be eliminated by the scaling property 
when a trapping potential vanishes. Section IlIII presents 
the main results of this paper, in which a BEC droplet 
is shown to be stabilized in 3D free space by use of oscil- 



lating interactions. The stability diagrams with respect 
to the interaction and dissipative parameters are also ob- 
tained. Section Hvl studies the variational analysis, and 
Sec. M concludes this paper. 



II. THE SCALING PROPERTY OF THE 
GROSS-PITAEVSKII EQUATION 

We consider the GP equation in the presence of dissi- 
pation given by |0, 0| 
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where 7 is a phenomenological dissipation constant which 
is to be determined experimentally, and the trapping po- 
tential Vtrap(r,i) and the s-wave scattering length a{t) 
are controlled to vary in time. The last term on the 
right-hand side of Eq. Q guarantees the normalization 
/ (ir|-0p = N with TV being the number of atoms, where 
the chemical potential is given by 



M = / drV/ -T^V" + ytrap(r, t) + 
' * 2m 



... 

For simplicity, we assume an isotropic trapping poten- 
tial 



Vtrap(r,t) — /trap(i)- 



(3) 



where /trap(0 is unity aX t — Q, then decreases towards 
zero, and vanishes at some time. The s-wave scattering 
length is made to oscillate as 



a{t) = /int(i)(ao + ai sinVLt), 



(4) 



where fint{t) ramps up from to 1. The ramp functions 
/trap(^) and /int(i) are introduced to avoid initial nona- 
diabatic disturbances that cause dynamic instabilities. 
These functions must be chosen so that neither collapse 
nor expansion occurs in the course of the ramping. 

Normalizing the length, time, energy, and wave func- 
tion by {h/mfl)^/^ = do, hfl, and ^/Nd^^^^, respec- 
tively, we obtain the normalized form of the GP equation: 
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where = ATTNa{t)/do = fint{t){go + gi sint), and the 
wave function is normaUzed as /drj'^p = 1. It should 
be noted that once the trapping potential is switched off 
[/trap(0 = 0], Eq. © no longer depends on U exphcitly. 
This implies that the CI dependence can be eliminated 
by appropriate scaling of the time, length, and interac- 
tion. The dynamics at driving frequency af2 can be re- 
duced to those at O by rescaling t — > at, r y^r, and 
a(t) — > ^/aa{t). For example, if large Q is unavailable 
experimentally, we have only to increase the strength of 
interaction and enlarge the size of the BEC droplet. 



III. NUMERICAL SIMULATIONS 

We numerically solve the GP equation ||SJ) using the 
Crank-Nicholson scheme ^3 : where the following simple 
forms of the ramp functions are used: 



/trap(i) 
/int(t) 



{t> Ttrap), 

1 - {t/Tint - 1)2 (0 < t < Tint) 
1 {t > Tint). 



(6) 
(7) 



We use a quadratic function for /int(i) because it makes 
initial nonadiabatic disturbances smaller than linear 
ramps. 

Figure n shows the time evolution of the peak density 
\ip{r = 0)p and the monopole moment (r) = J (irr|-0P 
for go = —69 and gi = 155. The dissipation constant 
7 is taken to be 0.03 [lolflll|. and the ramp parameters 
are Ttrap = 16 and Tint = 10. The initial state is chosen 
to be the noninteracting ground state for a trapping po- 
tential with uj = ri/30. The density and width oscillate 
at the driving frequency. [The black bands in Fig. ^ (a) 
represent rapid oscillations of (r) and |'(/'(0)P which are 
beyond the resolution limit.] We note that the ampli- 
tudes and mean values of both (r) and |V'(0)P converge 
to constant values, indicating that a BEC droplet is dy- 
namically stabilized. The initial disturbances last for a 
relatively long time shown as slow oscillations or ripples 
of (r) and |V'(0)P in Fig- Q (a), which gradually decay 
due to dissipation. From the time evolution of the den- 
sity profile |'0(?')P [inset in Fig. ^(b)], we find that the 
droplet exhibits breathing-mode oscillations with a large 
amplitude. This indicates that large oscillations of mass 
flow persist even after stabilization. We find from Fig. 
(b) that the phase of the breathing-mode oscillations of 
the droplet lags slightly behind that of the oscillating in- 
teraction shown as the dotted curve. This retardation is 
considered to be due to dissipation because there is no 
phase difference in the dissipation-free 2D case jsf , where 
|V'(0)|2 and (r) are in-phase and out-of-phase with the 
oscillating interaction. 
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FIG. 1: (a) Time evolution of the peak density |7/)(r — 0)f 
(left axis) and monopole moment (r) — J drr\tlj\^ (right axis) 
for the oscillating interaction g(t) — —69 -I- 155 sin t with 
7 = 0.3. The initial state is the noninteracting ground state 
in a trapping potential r^/1800. Then the interaction is grad- 
ually turned on, and the trapping potential is simultaneously 
turned off as shown in the inset, (b) A magnification of (a). 
The inset shows the density profile |i/)(r)|^ from t = 8000 to 
t — 8020. The dotted line is a sine function for comparison of 
the phase. 



We gradually changed go and gi from those of the sta- 
ble droplet state in Fig. ^ and found distinct types of 
instabilities as shown in Fig. |21 When \go\ is decreased, 
a droplet expands because of a decrease in an attractive 
interaction as shown in Fig. |21 (a) . The expansion occurs 
also by an increase in gi, since the effective repulsive in- 
teraction due to the oscillation is proportional to g^ [see 
Eq. ((TT]) ]. When \go\ is increased and gi is decreased, 
the system becomes dynamically unstable against the 
slow oscillation and eventually expands away as shown 
in Fig. 13 (b)- The expansion occurs, despite the fact 
that an increase in \gQ\ enhances the attraction and a 
decrease in gi suppresses the repulsion. The slow oscil- 
lation is seen in the profile of the maxima of [^(O)!^ (the 
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FIG. 2: Time evolution of the peak density \tp{r = 0)|^ (left 
axis) and monopole moment (r) = J drr\ip\^ (right ax;is). The 
initial state is prepared as in Fig. Q(a), and {go,gi) are lin- 
early changed from (-69, 155) to (-69, 160) in (a), (-69, 129) 
in (b), and (-72, 173) in (c) during <t < 600. The insets 
in (b) and (c) are magnifications of \ip{r = 0)p in the dashed 
squares and the density profiles [■(/'(r)^ at the times indicated 
by A and B. 



ripple in the upper edge), while the minima is almost 
constant. Figure |3 (c) shows an instability that arises 
for large \go\ and gi. Since the modulation at half fre- 
quency grows, this instability is considered to arise from 
parametric resonance. The density profiles in the inset 
show that higher radial modes are excited, which leads 
to destabilization of the droplet. 

We prepare a stable droplet in the same manner as 
in Fig. ^ and change the parameters to various values to 
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FIG. 3: The lifetime of a BEC droplet after a stable state is 
prepared as in Fig. Q and the interaction is linearly changed 
during < t < 600 as in Fig. |5| (a) The black region cor- 
responds to the one in which a droplet expands during the 
change of the interaction, and the white region corresponds 
to the one in which a droplet survives until at least t = 10600. 
The behaviors in the regions "instability I, II, and III" corre- 
spond to those in Fig.|5|(a), (b), and (c), respectively, (b) The 
parameter 7 is also linearly changed as well as the interaction 
during < t < 600 from 7 = 0.03 to 0.01. 



obtain the stability diagram shown in Fig.|31 The param- 
eters go and gi are changed linearly during < t < 600. 
We define the lifetime of the droplet as the duration be- 
tween t — 600 and the time at which the droplet begins to 
expand indefinitely. The regions referred to as "instabil- 
ity I, II, and III" in Fig. |3| correspond to the instabilities 
demonstrated in Fig. |21(a), (b), and (c). We find from 
Fig- ini {^) that the stable region is surrounded by the 
three distinct types of instabilities. The boundary be- 
tween the "stable" and "instability I" regions are sharp 
while the "instability 11" gradually sets in as \go \ increases 
and gi decreases. The collapse occurs for large \go\ and 
small gi, but the corresponding parameter regime lies 
outside Fig. O 

In Fig|31(b), 7 is linearly decreased from 0.03 to 0.01 
as \go\ and gi are (also linearly) changed. The instabil- 
ity regions enlarge as compared with those of Fig|21 (a), 
and the stability region disappears. From this result, we 
conclude that the stability region exists only for 7 7^ 0. 
However, we cannot exclude the possibility that other 
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FIG. 4: The profiles of the column density of a BEG droplet 
obtained by full 3D numerical calculations. Since the isotropic 
symmetry is preserved, we obtain the same column density 
from any direction. The size of the images is 12 x 12 in units 

of {h/mny^^. 



stable states exist, since we do not investigate the entire 
functional space. 

In numerical calculations of the time evolution, we 
must pay special attention to the boundary effect. Since 
the trapping potential is absent, the atoms that escape 
from the droplet spread out, reflect at the boundary, and 
return to the droplet region, producing spurious bound- 
ary effects. For example, when the Dirichlet boundary 
condition is imposed at r = 116, the lifetime of the 
droplet (whose definition is the same as in Fig. ^ for 
go = —69 and gi = 134 becomes ~ 4100, while the cor- 
rect value is ~ 3600. Therefore, the spatial cutoff in the 
numerical analysis must be much larger than the size of 
the droplet. 

The numerical simulations that we have presented so 
far have been carried out under the assumption of the 
spherical symmetry of the system. However, this ge- 
ometry cannot generate multipole dynamical instabili- 
ties which destroy a BEC droplet. In order to check 
the absence of such multipole instabilities, we have per- 
formed full 3D calculations by discretizing the space of 
100 X 100 X 100 size in our dimensionless unit into a 
256 X 256 X 256 mesh. A limited spatial size due to our 
computational power causes spurious boundary effects 
and a rough mesh produces numerical errors. Neverthe- 
less, as shown in Fig. 01 snapshots of time evolution of a 
BEC droplet appear isotropic and do not show any mul- 
tipole instability. We have thus confirmed that a BEC 
droplet can be stabilized without multipole instabilities. 



IV. VARIATIONAL ANALYSIS 

The Gaussian variational wave function well describes 
the dynamical stabilization of BECs in 2D free space ^ 
Q. We examine a variational method in 3D using the 
Gaussian trial function as [T^ll^ 

1 / r2 Rr^\ 

^va.(r, t) = -^j^ exp (8) 

where R{t) is the variational parameter that character- 
izes the size of the condensate. Substituting Eq. into 



the action 

S = J dtdvr - ^ + f IV-P) V', (9) 

we obtain the equation of motion for R as 
d f 1 G \ 

where G{t) = {go + gi smt)/{2^/'^TT^^^) = Go + Gi sint. 

We separate R into the slowly varying part Rq and the 
rapidly oscillating part p as R = Ro + p. According to 
Ref. 7], an effective potential for i?o is given by fi/{AQ'^), 
where /i is the amplitude of the oscillating "force" for Rq 
and is its frequency. From Eq. (|10|l . fi corresponds to 
Gi/(2_Rq). Since we are using a system of units in which 
57 = 1, the effective potential becomes G^/(16i?o); ^'^d 
then the equation of motion for Rq reads 

° dRo\2Rl^ 6Rl WR^o) ~ dRo ' ^ ' 

This equation agrees with that in Ref. i^j when fl is much 
larger than the characteristic frequencies of the system. 
Thus, the oscillating interaction blocks the collapse by 
the effective potential proportional to Rq^ ■ 

The effective potential Ucs has a local minimum for 
appropriate values of Go and Gi. For instance, when 
Go = -70/(2i/27r3/2) and Gi = 140/(21/^3/2), U,b has 
a local minimum at i?o — 2.4. In fact, by numerically 
solving the equation of motion (|10|) for R, we find that 
a stable solution exists. However, Fig.|3|(b) implies that 
the system is unstable against the dynamical instability 
for 7 = 0. This indicates that the simple Gaussian func- 
tion © is insufhcient to describe the instability of the 
system. 

Figure |31 illustrates the difference between the density 
profiles of the stable droplet {t > 8000 in Fig. ^ and the 
Gaussian functions fitted to them. We adopt the least 
squares fitting of r^\^p{r)\^ to /{■K^^'^a'^) exp{—r'^ /a'^) 
with fitting parameter a, which is an appropriate fitting 
in 3D. Figure El shows that the numerically exact wave 
functions significantly deviate from the Gaussian func- 
tions even in the stable droplet. We note that the Gaus- 
sian function in Fig |S1 (a) is almost the same as that in 
FigEl (b) . This indicates that the outer region, which has 
large weight in the fitting due to the factor r^, is almost 
stationary, and the mass flow occurs mainly around the 
center of the droplet. Thus, more appropriate functions 
are needed to correctly describe the phenomena. 



V. CONCLUSIONS 

We have studied a Bose-Einstein condensate with an 
oscillating interaction with dissipation, and shown that a 
BEC droplet is stabilized in 3D free space with realistic 
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FIG. 5: The density profiles (solid curves) of the 

stable droplet for go = -69 and gx = 155 (t > 8000 in Fig.IU 
when the peak density becomes (a) maximal and (b) minimal. 
The dashed curves are the Gaussian functions fitted to the 
density profiles. 

dissipation, 7 = 0.03. The frequency of the osciUating in- 
teraction must be much faster than the characteristic fre- 
quency of the system, but it can be a moderate frequency 
according to the scaling property discussed in Sec. 
For instance, the situation in Fig. ^ can be realized in a 
condensate of 10'* ®^Rb atoms by the s-wave scattering 
length a{t) = -0.6 + 1.35 sin f7i [nm] with = 100 x 27r 
Hz, where a trapping potential with frequency 3.3 Hz 
is removed in 25 ms, and the final size of the droplet 
becomes a few micrometers. Such an oscillation of in- 



teraction can be easily realized in experiments using the 
Feshbach resonance |l5| . 

Under gravity, a condensate falls after the trapping 
potential is turned off. In the above example, the grav- 
itational sag in the initial trapping potential is ~ 2 cm, 
and the condensate falls about 1 mm until the trapping 
potential vanishes. The effect of gravity can be canceled 
out using the magnetic levitation 16], which enables us 
to observe the longtime behavior of the droplet. 

A "gaseous EEC droplet" found in this paper is dif- 
ferent from the usual condensate in that it coheres by 
itself without the help of the trapping potential. This 
self-trapped matter wave might exhibit interesting dy- 
namics, such as in collective mode, collapsing dynamics, 
and vortex nucleation, which will be discussed elsewhere. 

Note added — While preparing the present paper, a 
similar paper appeared which also Reports stabiliza- 
tion of a BEC in 3D free space by an oscillating interac- 
tion. The strength of interaction studied in Ref. 0| is 
much larger than that of ours. In our parameter regime, 
a droplet can be stabilized adiabatically by simple ramp 
functions with realistic dissipation. On the other hand, in 
Ref a droplet is stabilized without dissipation by fol- 
lowing a complicated ramp scheme to reach a final state. 
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